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Neutron stars—compact objects with masses similar to that of our Sun but radii comparable to the size of a 
city—contain the densest form of matter in the universe that can be probed in terrestrial laboratories as well as in 
earth- and space-based observatories. The historical detection of gravitational waves from a binary neutron star 
merger has opened the brand new era of multimessenger astronomy and has propelled neutron stars to the center 
of a variety of disciplines, such as astrophysics, general relativity, nuclear physics, and particle physics. The 
main input required to study the structure of neutron stars is the pressure support generated by its constituents 
against gravitational collapse. These include neutrons, protons, electrons, and perhaps even more exotic con- 
stituents. As such, nuclear physics plays a prominent role in elucidating the fascinating structure, dynamics, and 


composition of neutron stars. 


I. INTRODUCTION 


Among the most interesting questions animating nuclear 
astrophysics today are: What are the new states of matter at 
exceedingly high density and temperature? and how were the 
elements from iron to uranium made? In one clean sweep, 
the historical detection of gravitational waves from the binary 
neutron star merger GW 170817 by the LIGO- Virgo collabora- 
tion [LL] has provided critical insights into the nature of dense 
matter [2/19] and on the synthesis of the heavy elements in 
the cosmos [2023]. 

The term neutron star appears in writing for the first time 
in the 1933 proceedings of the American Physical Society by 
Baade and Zwicky who wrote: With all reserve we advance 
the view that supernovae represent the transition from ordi- 
nary Stars into “neutron stars”, which in their final stages 
consist of extremely closed packed neutrons (24). It appears, 
however, that a couple of years earlier Landau speculated on 
the existence of dense stars that look like giant atomic nu- 
clei [25]. The first quantitative calculation of the structure 
of neutron stars was performed by Oppenheimer and Volkoff 
in 1939 by employing the full power of Einstein’s theory of 
general relativity [26]. Using what it is now commonly re- 
ferred to as the Tolman-Volkoff-Oppenheimer (TOV) equa- 
tions [26}|27], Oppenheimer and Volkoff predicted that a neu- 
tron star supported exclusively by quantum mechanical pres- 
sure from its constituent neutrons will collapse into a black 
hole once its mass exceeds seven tenths of a solar mass, or 
Mmax =0.7Me. Almost 30 years later, Jocelyn Bell discov- 
ered “pulsars” which, after a period of great confusion in 
which they were mistaken as potential beacons from an ex- 
traterrestrial civilization, were finally identified as rapidly ro- 
tating neutron stars [28]. 

Since then, the field has evolved by leaps and bounds, 
and to date, a few thousands neutron stars have been identi- 
fied [29]. Among the most impactful discoveries since then 
is the identification of three neutron stars with masses in the 
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vicinity of two solar masses [30}{33}, although a recent publi- 
cation has reported the observation of a neutron star with an 
even larger mass of about Mmax =2.35Mo [34]. These dis- 
coveries suggest that the original prediction by Oppenheimer 
and Volkoff underestimates the current limit on the maximum 
neutron star mass by at least a factor of three. It is now known 
that such dramatic discrepancy is associated with the unreal- 
istic assumption that neutrons behave as a non-interacting gas 
of fermions. That is, Oppenheimer and Volkoff incorrectly as- 
sumed that neutrons in a neutron star behave as electrons in a 
white-dwarf star. Thus, the mere existence of massive neutron 
stars highlights the vital role of nuclear interactions in explain- 
ing the structure, dynamics, and composition of neutron stars. 
In this manner, the detailed understanding of neutron stars has 
been transferred to the domain of nuclear physics. 

The TOV equations are model independent as the only as- 
sumption behind them is the validity of Einstein’s theory of 
general relativity. However, the solution of the TOV equa- 
tion requires a detailed description of the relation between the 
energy density—which provides the source of gravity—and 
the pressure—which prevents the gravitational collapse of the 
star. This relation is enshrined in the equation of state (EOS). 
Unlike a classical ideal gas where the density is low and the 
temperature is high, the relevant domain for neutron stars in- 
volves high densities and low temperatures. That is, in neu- 
tron stars the inter-particle separation is small relative to the 
thermal de Broglie wavelength, thereby highlighting the role 
of quantum correlations. That the inter-particle separation in 
the core of neutron stars is small implies that the dynamics 
is dominated by the short-distance repulsion of the nucleon- 
nucleon interaction. Such a repulsion provides an additional 
source of pressure support above and beyond the Fermi (or 
Pauli) pressure assumed by Oppenheimer and Volkoff. It is 
this additional support that becomes essential in supporting 
massive neutron stars against gravitational collapse. Further, 
the quantum regime of low temperatures and high densities 
suggests that neutron stars may effectively be treated as zero- 
temperature systems. Thus the challenge is to infer the ab- 
solute ground state of the system over the enormous range of 
densities and pressures present in a neutron star. 


Il. TOMOGRAPHY OF A NEUTRON STAR 


The structure of neutron stars is both interesting and com- 
plex. Figure{I] encapsulates our current understanding of the 
structure and composition of neutron stars [354/37]. The out- 
ermost surface of the neutron star contains a very thin atmo- 
sphere of only a few centimeters thick that is composed of Hy- 
drogen, but may also contain heavier elements such as Helium 
and Carbon. The electromagnetic radiation emitted from the 
surface and collected on earth- and spaced-based telescopes 
is instrumental in our determination of some critical param- 
eters of the neutron star. For example, under the assumption 
of blackbody emission from the stellar surface, one can deter- 
mine the effective surface temperature T of the star. More- 
over, by invoking the Stefan-Boltzmann law, L = 4noR?T*, 
one can in principle infer the radius of the star, a quantity that 
provides critical information on the equation of state. Un- 
fortunately, complications associated with distance measure- 
ments that are required to determine the absolute stellar lumi- 
nosity L, as well as atmospheric distortions of the black-body 
spectrum, make the accurate determination of stellar radii a 


challenging task [38}/40)}. 
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Figure 1. A theoretically-informed rendition of the structure and 
phases of a neutron star. Courtesy of Dany Page [36]. 


Fortunately, the situation has improved through a better 
understanding of systematic uncertainties, important theoreti- 
cal developments, and the implementation of robust statistical 
methods [41}/48]. More recently, the tidal deformability of 
neutron stars extracted from GW170817 together with 
the monitoring of hot spots in neutron stars [49153], are pro- 
viding new tools for the accurate determination of stellar radii. 

Just below the atmosphere lies a thick envelope that acts 
as a “blanket” between the hot interior and the relatively cold 
surface [35]. Matter in this region is not yet fully degener- 
ate so details associated with both the complexity of the at- 
mosphere and large magnetic fields are beyond the nuclear 
physics purview. However, the layers below, starting with the 


solid crust and ending with the high density core, are the main 
focus of this article. The non-uniform crust represents a region 
of about 1 km that is composed of a lattice of exotic neutron- 
rich nuclei immersed in a uniform, neutralizing electron gas 
and a quantum liquid of superfluid neutrons [55]. Even 
deeper in the crust, the universal phenomenon of Coulomb 
frustration favors the formation of a complex set of topologi- 
cal phases dubbed “nuclear pasta’ [56] [57]. The non-uniform 
crust sits above a uniform liquid core that consists of neutrons, 
protons, electrons, and muons. The core accounts for practi- 
cally all the mass and for about 90% of the size of a neutron 
star. Finally, depending on the presently unknown largest den- 
sity that may be attained in the stellar core, there is also a pos- 
sibility, marked with a question mark in Fig{I| of a transition 
into a phase made of exotic particles, such as hyperons, meson 
condensates, and deconfined quark matter [36]. 


A. The Outer Crust 


The typical densities encountered in the outer crust of the 
neutron star span about seven orders of magnitude, from 
10*g/cm? to about 10!'g/cm? [58]. In this region the average 
inter-particle separation is significantly larger than the range 
of the nucleon-nucleon interaction. As such, the system—if 
arranged in a uniform configuration—would not be able to 
benefit from the intermediate-range attraction of the nucleon- 
nucleon interaction. Hence, at these densities it becomes ener- 
getically favorable for the system to break the uniform spatial 
symmetry and to cluster into neutron-rich nuclei segregated 
into a Coulomb lattice that is immersed in a uniform electron 
gas. Electrons are an essential component of the neutron star 
as they enforce the overall charge neutrality of the system. 
Thus, the constituents of the outer crust are neutrons, protons, 
and electrons in chemical (or beta) equilibrium. 

In the outer crust, and indeed throughout the entire star, the 
basic question that one needs to answer is how do the relevant 
constituents organize themselves to minimize the energy over 
the enormous range of densities encountered in a neutron star. 
Fortunately, the dynamics of the outer crust is relatively sim- 
ple, as its dependence on the strong nuclear interaction is fully 
encapsulated in the masses of a few atomic nuclei. The dy- 
namics of the outer crust is entirely contained in three distinct 
contributions to the Gibbs free energy of the system: lattice, 
electronic, and nuclear. The energy associated with the lat- 
tice contribution is complicated due to the long-range nature 
of the Coulomb interaction. It consists of divergent contribu- 
tions that must be carefully canceled as required by the overall 
charge neutrality of the system. However, accurate numerical 
calculations for the electron gas have been available for a long 
time [59] [60], which have since been adapted to the problem 
of interest [62]. In turn, the electronic contribution 
is encoded in the known energy per particle of a relativistic 
Fermi gas of electrons. All that is left is the nuclear contri- 
bution. The nuclear contribution is obtained by identifying 
the optimal nucleus whose mass per nucleon minimizes the 
Gibbs free energy per particle at a given pressure. In essence, 
all one needs to compute the composition of the outer crust is 


an accurate nuclear mass table. 
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Figure 2. A theoretically-informed rendition of a neutron star with 
emphasis on the crustal composition. Courtesy of Sanjay Reddy. 


At the top layers of the outer crust where the density is low- 
est, it is energetically favorable for nucleons to cluster into 
56Re nuclei, which arranged themselves in a face-centered cu- 
bic lattice [58]; see Fig..2| Iron is the optimal nucleus at low 
density because it has the lowest mass per nucleon in the entire 
nuclear chart. However, as the density and pressure increase, 
Fe ceases to be the preferred nucleus. The reason for this 
behavior is that the electronic contribution to the Gibbs free 
energy grows rapidly with increasing density/pressure. Thus, 
it becomes energetically advantageous for the energetic elec- 
trons to capture onto protons, making the matter neutron rich. 
At densities of about 10°g/cm?, °°Ni becomes the most fa- 
vorable nucleus. Note that although °°Ni has a lower bind- 
ing energy per nucleon than *Fe, it is °°Fe which has a lower 
mass per nucleon because its slightly smaller neutron fraction; 
the entire difference amounts to only 12 keV. For a neutron 
star in hydrostatic equilibrium, the pressure and density con- 
tinue to increase as one moves towards the deeper layers of 
the outer crust. As the density increases, the system evolves 
into a Coulomb lattice of progressively more exotic, neutron- 
rich nuclei [61] (62). Finally, at a neutron-drip density of about 
10!!g/cm3, atomic nuclei are incapable of binding any more 
neutrons. As suggested in Fig..2| most of the existing mass 
models predict that the sequence of progressively more ex- 
otic nuclei ends with ''’Kr—a nucleus with 36 protons and 
82 neutrons. 

Ultimately, the determination of the optimal nucleus is 
set by a competition between the electron fraction and the 
nuclear symmetry energy, a concept that appears for the 
first time in the semi-empirical mass model of Bethe and 
Weizsicker [63] [64]. The symmetry energy quantifies the en- 
ergy cost of turning symmetric nuclear matter, with equal 
number of protons and neutrons, into pure neutron matter. 
Low electron fractions diminish the electronic contribution 
to the Gibbs free energy but at the expense of creating very 
neutron-rich nuclei. In turn, a large neutron excess is dis- 
favored by the nuclear symmetry energy. Knowledge of the 
symmetry energy—and especially its density dependence—is 
one of the most important contributions of nuclear physics to 


our understanding of neutron stars. 

Given that the last krypton isotope with a well measured 
mass is °’Kr, which is still 21 neutrons away from the pre- 
dicted drip-line nucleus !!8Kr, one must rely on mass models 
that are often hindered by uncontrolled extrapolations. To mit- 
igate this problem, two important developments have taken 
place. First, the highly anticipated Facility for Rare Isotope 
Beams (FRIB) is now up and running. FRIB will surpass cur- 
rent limits of precision and sensitivity that will enable the pro- 
duction of exotic nuclei at or near the drip lines. Second, state- 
of-the-art mass models have been refined by capital- 
izing on the explosion of machine learning techniques 
[82]. To overcome some of the intrinsic limitations of existing 
mass models, machine learning approaches have been imple- 
mented with the goal of optimizing mass residuals between 
theory and experiment. Beyond a significant improvement in 
the predictions of existing models, the Bayesian refinement 
provides predictions that are accompanied by statistical un- 
certainties. 

Finally, we note that there are three important regions in 
the nuclear chart that are critical for our understanding of the 
outer crust: (a) the iron-nickel region, (b) the region around 
the N =50 isotones, and (c) the region around the N = 82 iso- 
tones. In the case of the iron-nickel region, the masses of 
all relevant nuclei have been measured with enormous pre- 
cision. Although not with the same level of precision as in 
the iron-nickel region, nuclear masses in the neighborhood 
of the N = 50 isotonic chain have for the most part been 
measured [83]. For example, in a landmark experiment at 
ISOLTRAP at CERN, the mass of **Zn (Z =30 and N =52) 
was measured for the very first time [84]. The addition of this 
one mass value alone resulted in an interesting modification to 
the composition of the outer crust [85]. Lastly, very little 
(if at all!) is known about the masses of the very neutron-rich 
nuclei in the neighborhood of the N = 82 isotonic chain. It 
is precisely in this region that machine learning refinements 
to mass models informed by experiments at FRIB and other 
rare-isotope facilities will become of outmost importance. 


B. The Inner Crust 


The inner stellar crust comprises the region from neutron- 
drip density up to the density at which uniformity in the sys- 
tem is restored. The exact crust-to-core transition density 
is unknown, but calculations predict that it lies between a 
third to a half of nuclear matter saturation density. Based 
on the nearly uniform interior density of heavy nuclei, nu- 
clear matter saturation density has been determined to be 
around n, =0.15 fm~? [86], corresponding to a mass density 
of p, + 10!4g/cm?. 

On the top layers of the inner crust, nucleons continue to 
cluster into a Coulomb crystal of neutron-rich nuclei embed- 
ded in a uniform electron gas. However, beyond neutron drip, 
the crystal is in chemical equilibrium with a superfluid neu- 
tron vapor. Most of the evidence in favor of neutron pairing in 
the crust is obtained from studying its early thermal relaxation 
period [87]. An interesting consequence of superfluid pairing 


is the appearance of pulsar “glitches”. Rotation-powered pul- 
sars tend to spin down slowly and steadily due to the emission 
of magnetic dipole radiation, making pulsars one of nature’s 
most accurate clocks. For example, the 33 ms Crab pulsar 
spins down by about 13 us per year. Despite this remarkable 
regularity, young pulsars often glitch, a unique phenomenon 
characterized by a sudden and abrupt spin-up in their rota- 
tional frequency. Pulsar timing has revealed that glitches are 
recurrent, with some of the more active glitchers being the 
Vela and the Crab pulsars [87]. Although the precise details of 
the pulsar glitch mechanism are unclear [88] [89], it is believed 
that the glitch mechanism is intimately related to the forma- 
tion of superfluid vortices in the inner crust of the rotating neu- 
tron star; see Refs.(90H92] and references contained therein. 
As the pulsar slows down, the initial distribution of vortices 
that are believed to be pinned to the crystal lattice falls out of 
equilibrium. This induces a differential rotation between the 
slower neutron star and the faster superfluid vortices. When 
the differential lag is too large, then suddenly and abruptly 
some fraction of the vortices unpin, migrate outwards, and 
transfer their angular momentum to the solid crust—which is 
detected as a glitch. After this transient event, the star contin- 
ues to slow down, stresses between the crust and the superfluid 
start to build up again until eventually more vortices unpin, 
transfer angular momentum to the solid crust, and ultimately 
generate another glitch. 


From the perspective of nuclear structure, the most fasci- 
nating component of the inner crust is the nuclear pasta phase. 
In the outer crust and in the top layers of the inner crust, the 
short-distance scales associated with the strong interaction is 
well separated from the large distance scales characteristic 
of the Coulomb repulsion. In this density regime, nucleons 
bind onto nuclei that are well segregated in a crystal lattice. 
However, at the bottom layers of the inner crust, these length 
scales become comparable. Competition between short-range 
attraction and long-range repulsion leads to a universal phe- 
nomenon known as “Coulomb frustration”. Such a competi- 
tion becomes responsible for the emergence of complex struc- 
tures of various topologies that are collectively known as nu- 


clear pasta [56] [57]. 


Most theoretical approaches to the structure and dynamics 
of the inner crust fall into two broad categories: (a) Mean-field 
models using accurately calibrated interactions that incorpo- 
rate quantum effects but miss some important many-body cor- 
relations [944100] and (b) semiclassical Monte Carlo (MC) 
and Molecular Dynamics (MD) formulations where many- 
body correlations are treated properly but quantum effects are 
incorporated in an approximate way [93] [101//108]. Among 
the virtue of the semi-classical MC/MD approaches is that ex- 
isting computational resources afford the possibility of simu- 
lating an enormous number of particles, thereby minimizing 
the finite-size effects. For example, Schneider and collabora- 
tors have studied the formation of multi domains in the nu- 
clear pasta using a simulation box containing about three mil- 
lion nucleons [106]. Figure[3]displays a snapshot of an early 
Monte-Carlo simulation for a system of 4,000 nucleons at a 
density that clearly identifies the formation of non-spherical 
clusters embedded in a dilute neutron vapor [93] [109]. An- 


Figure 3. A snapshot of a Monte-Carlo simulation for a system con- 
taining 4000 nucleons (3200 neutrons—800 protons) at a baryon den- 
sity of n=n,/6 and an effective temperature of T=1 MeV , 


other great virtue of the semiclassical simulations is that pasta 
formation is studied in an unbiased way without assuming any 
particular set of shapes. Rather, the system evolves dynami- 
cally into these complex configurations from an underlying 
two-body interaction consisting of a short-range nuclear at- 
traction and a long-range Coulomb repulsion. This fact is 
of great relevance given that one of the main characteristics 
of Coulomb-frustrated systems, such as the nuclear pasta, is 
the preponderance of quasi-degenerate ground states. That 
is, small changes in either the baryon density or the proton 
fraction can induce dramatic changes in topology, but with 
only a very modest gain in energy [110]. However, such small 
changes may have a profound impact on the crustal dynamics 
and transport properties [93] [109] [TI 1/119]. 

Finally, whereas the study of nuclear pasta is fascinating in 
its own right, it may be the key to explain some interesting as- 
trophysical phenomena such as the lack of X-ray-emitting iso- 
lated pulsars with long spin periods and the slow cooling 
rate of MXB 1659-29 [121], a binary system consisting of a 
neutron star accreting from a low-mass companion star. How- 
ever, it is unclear how thick the crustal region is that is occu- 
pied by the nuclear pasta. This feature—as well as the transi- 
tion density from the non-uniform crust to the uniform core— 
depend critically on the poorly known density dependence of 
the symmetry energy [122]. If the symmetry energy is small 
at the densities of relevance for pasta formation, then the im- 
balance between neutrons and protons may be large, thereby 
reducing the overall long-range repulsion and thus hindering 
the emergence of the pasta phase; after all, a system made 
entirely of neutrons does not cluster. 


C. The Outer Core 


At densities of about 10!'*g/cm? the nuclear pasta “melts” 
and uniformity in the system is restored. Above the crust- 


core transition density the original perception of Baade and 
Zwicky of a neutron star is finally realized: a uniform 
assembly of extremely closed packed neutrons. At these den- 
sities the neutron star consists of neutrons, protons, electrons, 
and muons in chemical equilibrium. As a result of chemical 
equilibrium, the proton fraction in the core hovers around 10% 
with a corresponding fraction of charged leptons required to 
maintain charge neutrality. Muons are important once the 
Fermi energy of the electron equals the rest mass of the muon, 
which happens in the vicinity of nuclear matter saturation den- 
sity. 


Whereas the stellar crust displays fascinating and intriguing 
dynamics, its structural impact on the star is rather modest. In- 
deed, with a crustal radius of about 1-2 km, the core accounts 
for about 90% of the size and most of the mass of the neutron 
star. Thus, the stellar core provides a unique laboratory for the 
study of the ground state of uniform, neutron-rich matter over 
conditions that are unattainable in terrestrial laboratories. By 
the same reason, the equation of state of neutron-rich matter is 
poorly constrained by laboratory observables. As articulated 
later in Sec{III] the most powerful constraints on the EOS be- 
tween the core-crust transition density and about twice satu- 
ration density emerge from a combination of theory and lab- 
oratory experiments. The stellar composition in this region 
is also controlled by the density dependence of the symmetry 
energy. A stiff symmetry energy, namely, one that increases 
rapidly with baryon density, disfavors a large neutron-proton 
asymmetry. This has interesting consequences as a stiff sym- 
metry energy favors large proton fractions that are required 
for the onset of the direct Urca process [123]. The direct Urca 
process provides a very efficient mechanism for cooling the 
star by neutrino emission from direct beta decay and electron 
capture [351/37] [124]. A recent laboratory experiment that ex- 
tracted the neutron skin thickness of 7°*Pb seems to suggest 
that the symmetry energy in the vicinity of nuclear satura- 
tion density is stiff [125]. One of the many implications of 
this experiment is the possibility that the onset of direct Urca 
process happens at relatively low central densities and corre- 
spondingly for relatively low mass neutron stars [126]. This 
is particularly interesting and timely given that x-ray observa- 
tions suggest that some neutron stars—for example, the neu- 
tron star in the low-mass x-ray binary MXB 1659-29—may 
need a fast neutrino-cooling process [121]. Although there is 
still much to learn about the dynamics of neutron rich mat- 
ter in terms of conventional degrees of freedom, perhaps even 
more exciting is the possible emergence of exotic states of 
matter, besides the nuclear pasta and Coulomb crystals that 
populate the stellar crust. At densities above twice satura- 
tion density, a region denoted by a question mark in Fig{1| 
it is plausible that a system of nucleons and charged leptons 
evolves into others form of matter, such as hyperonic matter, 
meson condensates, deconfined quark mater, and/or color su- 
perconductors. In some cases, the argument in favor for such 
a transition is based on simple qualitative arguments but in 
other cases it requires detailed calculations. 


D. The Inner Core 


Qualitatively, the presence of hyperons in the core of neu- 
tron star is easy to understand. Hyperons are baryons contain- 
ing strange quarks and are therefore heavier than the neutrons 
and protons that populate the outer core. About 20% more 
massive than a neutron, the A-hyperon is a neutral baryon with 
one each of an up, a down, and a strange quark. At densities 
exceeding those found in the outer core, it becomes energeti- 
cally favorable to produce A-hyperons rather than to continue 
to add neutrons to the system. The emergence of hyperons, 
analogous to the emergence of muons, is thus a simple con- 
sequence of the Pauli exclusion principle. Naively, the onset 
of hyperons should happen around six times nuclear-matter 
saturation density. However, unlike leptons in high density 
matter, both nucleons and hyperons are subject to strong in- 
teractions in the stellar medium, which are poorly constrained 
at the short distances of relevance to the inner core. Thus, 
whereas it is plausible to expect the formation of hyperonic 
matter in the inner stellar core, the quantitative details of the 
transition are highly uncertain. Given that transitions to a 
new state of matter are usually accompanied by a softening of 
the equation of state, the emergence of hyperons yields max- 
imum neutron-star masses that are inconsistent with observa- 
tion [30133]. This problem is commonly referred to in the 
literature as the “hyperon puzzle” 127/137]. 


Among the other exotic states of matter that have been pro- 
posed to populate the inner stellar core are meson conden- 
sates, particularly pion and kaon condensates 
147]. Pions, kaons, and in general mesons are bound states 
containing one quark and one antiquark. The explanation for 
the possible emergence of meson condensates in neutron stars 
is relatively simple. The propagation of a pion or a kaon in the 
stellar medium can be drastically different from its propaga- 
tion in free space because of its coupling to “particle-hole” ex- 
citations. Unlike most atomic transitions, nuclear excitations 
often lead to collective modes involving the coherent superpo- 
sition of many particle-hole pairs. If the effective particle-hole 
interaction in the medium becomes strong and attractive, the 
assumed ground state may become unstable against the for- 
mation of pion- or kaon-like excitations. Such an instability 
signals the restructuring of the ground state into a new state 
of matter in which a fundamental symmetry of the underly- 
ing theory is no longer realized in the ground state. In this 
particular case, the broken symmetry is parity and the new 
ground state is a pion or a kaon condensate. Given the dif- 
ficulty in identifying robust astronomical signatures of such 
exotic states of matter and inspired by early theoretical calcu- 
lations of pion-like excitations in atomic nuclei [148}{150], ex- 
periments have searched—unsuccessfully—for precursors of 
pion condensation in laboratory experiments [1514154]. As 
in the case of hyperons, the difficulty in identifying the emer- 
gence of such exotic states of matter is our poor knowledge of 
the strong interaction under the extreme conditions present in 
the stellar core. 

Lastly, we address the possibility of matter composed of 
quarks and gluons populating the inner stellar core. Quarks 
and gluons are the basic constituents of Quantum Chromo Dy- 


namics (QCD) the fundamental theory of the the strong inter- 
action. On purely geometrical grounds one would expect a 
transition from hadronic matter, in which the quarks and glu- 
ons are confined within baryons and mesons, to deconfined 
quark matter. Given that protons and neutrons have a finite 
size, one would anticipate a transition to some form of quark 
matter around a critical density (7,) at which the nucleons start 
to touch. Assuming a close packing of spheres of radius 70, 
one obtains a surprisingly low critical density of 
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where the expression in the second bracket (~ 0.741) was 
shown by Gauss in 1831 to be the highest packing fraction, 
and as an estimate of 79 we used the charge radius of the 
proton ro + 0.84fm. At 279, QCD is notoriously difficult 
to solve in the non-perturbative regime, even using some of 
the most sophisticated theoretical approaches available today. 
Perhaps a more useful length scale is that at which the av- 
erage inter-quark separation becomes smaller than the QCD 
confining scale, so quarks could roam freely throughout the 
stellar core [29]. However, the naive picture of free-roaming 
quarks must be revisited because of the attractive interac- 
tion between a pair of quarks [155]. As shown by Bardeen, 
Cooper, and Schrieffer in their seminal work on superconduc- 
tivity, an arbitrarily weak attractive force may cause a dra- 
matic restructuring of the normal ground state. Under such 
a paradigm, QCD predicts that at ultra-high densities and 
low temperatures—where the up, down, and strange quarks 
are effectively massless—the ground state is a color super- 
conductor with a unique “color-flavor-locking” (CFL) pair- 
ing scheme [1551159]. In the CFL phase the exact color 
symmetry and the approximate flavor symmetry become in- 
tertwined [155]. Moreover, in this limit of massless quarks, 
three-flavor quark matter is automatically neutral without the 
need to add neutralizing leptons; indeed, “no electrons are re- 
quired and none are admitted” [158]. Unfortunately, it is now 
believed that the extreme densities required for the CFL phase 
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to emerge can not be reached in the stellar cores. So assess- 
ing the impact of QCD at the densities of relevance to neutron 
stars remains an important challenge that is starting to be ad- 
dressed in the regime in the very high-density regime where 
QCD is amenable to perturbation theory [160]. 


Il. TAMING GRAVITY: THE EQUATION OF STATE 


The structure of neutron stars emerges from a fierce com- 
petition between gravity and the internal pressure support that 
prevents the collapse of the star. In hydrostatic equilibrium, 
these two opposing forces are perfectly balanced. The law 
of universal gravitation formulated by Newton in 1687 posits 
that two massive particles attract each other with a force that 
is proportional to the product of their masses and inversely 
proportional to the square of their separation. From all fun- 
damental forces of nature—gravity, electromagnetism, weak, 
and strong—only gravity is always attractive. Such a picture 
finds its natural explanation in Einstein’s theory of general rel- 
ativity that attributes the attraction to the curvature of space- 
time that all objects—both massive and massless—experience 
as they travel through space. Although Newton’s law of uni- 
versal gravitation is adequate for the description of the mo- 
tion of planets around the Sun, the gravitational effects around 
a neutron star are so strong that modifications to Newtonian 
gravity are required. 


A. The Tolman-Oppenheimer-Volkoff Equations 


The structure of neutron stars is encoded in the Tolman- 
Volkoff-Oppenheimer (TOV) equations [26] |27], which repre- 
sent the generalization of Newtonian gravity to the domain of 
general relativity. For static, spherically symmetric stars in 
hydrostatic equilibrium the TOV equations may be written as 
a pair of coupled, first order differential equations: 


dP(r) a : ene . (2a) 
dr c PI —2GM(r) /c?r) 7 

ES) say, 28(r) > 4nr*p(r), (2b) 
dr 


where G is Newton’s gravitational constant and c is the speed 
of light in vacuum. The structural information is contained 
in M(r), e(r), P(r), and p(r), which represent the enclosed 
mass, energy density, pressure, and mass density profiles, re- 
spectively. The arrows in the above expressions represent the 
Newtonian limit of the TOV equations. These limits are ob- 


tained by assuming that the corrections from general relativity 
are negligible, as in the case of white-dwarf stars [161/164]. 


Once two boundary conditions are provided in terms of a 
central pressure P(0) =P. and an enclosed mass at the origin 
M(0)=0, the TOV equations may be solved numerically us- 
ing a suitable differential-equation solver, such as the Runge- 


Kutta algorithm [165]. The solution of the TOV equations 
over a wide range of central pressures can then be used to 
generate the entire mass-radius relation, often regarded as the 
“holy grail”’of neutron star structure. Once these equations 
are solved, the radius R and mass M of the neutron star are 
obtained as P(R) =0 and M=M{(R). That is, the radius is 
determined as the point at which the pressure goes to zero 
and the total mass as the value of the enclosed mass at such 
radius. We note in passing that with the exception of a few 
(mostly academic choices), the TOV equations must be solved 
numerically. This, however, presents some non-trivial chal- 
lenges associated with the diversity of scales involved in the 
problem. Neutrons, with a mass of m2 x 10~*’kg, provide 
most of the pressure support against the gravitational collapse 
of a neutron star with a typical mass comparable to that of 
our own Sun Mp ~2x10*9kg. This represents a mismatch 
in mass scales of 57 orders of magnitude. Thus, in order to 
avoid computer overflows—and as important, to gain critical 
insights into the problem—one must properly rescale the TOV 
equations [164]. It is interesting to point out that there are ap- 
proximately 10°’ neutrons in a neutron star, an astronomical 
number that makes Avogadro’s number of ~ 107? pale in com- 
parison! 

The validity of general relativity and the existence of a 
source of pressure support is all that is needed to derive the 
TOV equations. However, to solve the equations one must in 
addition provide an equation of state to connect the pressure to 
the energy density at zero temperature. The best known EOS 
is that of a classical ideal gas. However, unlike the low densi- 
ties and high temperatures that define a classical ideal gas, 
neutron stars are highly-degenerate compact objects where 
quantum effects are critical. In the quantum regime appli- 
cable to neutron stars, the average inter-particle separation is 
considerably smaller than the thermal de Broglie wavelength, 
so the fermionic (and bosonic if relevant) nature of the con- 
stituents plays a critical role. Moreover, whereas treating the 
constituents as a non-interacting Fermi gas may be adequate 
in the description of white-dwarf stars, assuming that neu- 
trons/protons in the stellar interior may be described as a non- 
interacting Fermi gas is woefully inadequate. Indeed, it was 
precisely such an assumption by Oppenheimer and Volkoff 
that led them to predict a maximum neutron star mass of only 
Mmax =0.7Mo [27], a prediction that is known to be in 
stark disagreement with the observation of = 2M. neutron 
stars [30433]. This underscores the critical role that strong 
interactions play in the description of neutron stars. In what 
follows, we discuss the various sources of pressure support 
and how recent theoretical, experimental, and observational 
advances are providing a new windows into the EOS. 


B. The Equation of State Ladder 


Given that neutron stars would inevitably collapse into 
black holes in the absence of pressure support, a detailed un- 
derstanding of the equation of state is critical. To do so, it 
is convenient to frame the discussion in terms of an “EOS 
ladder” akin to the cosmic distance ladder used in cosmol- 


ogy. In Fig [4] each rung in the ladder represents a theoreti- 
cal, experimental, or observational technique that informs the 
EOS in a suitable density regime. The first rung in the lad- 
der invokes chiral effective field theory, a purely theoretical 
approach rooted in the symmetries of QCD and which 
has its strongest impact on the EOS at and below saturation 
density. The next rung in the ladder includes laboratory ex- 
periments that constrain the EOS in the vicinity of saturation 
density. Beyond twice saturation density but still below the 
maximum density found in neutron stars, the EOS is informed 
by both electromagnetic observations and gravitational wave 
detections. For the highest densities anticipated to exist in the 
stellar core, the most stringent constraints are obtained from 
the identification of the most massive neutron stars. Finally, 
the last rung in the ladder is motivated by a recent claim that 
perturbative QCD may impose some constraints on the EOS 
of neutron-star matter [160]. The result is somehow surpris- 
ing given that QCD is estimated to become perturbative above 
40 times nuclear matter saturation density. Still, perturbative 
QCD may be effective in excluding some equations of state 
displaying some extreme behavior. In the next few sections 
we elaborate on the role of each each rung in the ladder. First, 
however, we start by explaining the marginal role that nuclear 
interactions play in generating the pressure support in the neu- 
tron star crust. 
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Figure 4. The nuclear equation of state ladder. Akin to the cos- 
mic distance ladder in cosmology, the equation of state ladder rep- 
resents a succession of theoretical, experimental and observational 
methods to determine how the pressure changes with increasing den- 
sity. Given that the range of densities spanned in a neutron star is 
enormous, no single method can determine the entire EOS. Instead, 
each rung of the ladder provides information that can be used to de- 
termine the EOS at a neighboring rung. 


Although the structure of the crust with all its exotic phases 
is fascinating, the role of nuclear degrees of freedom and their 
interactions on supporting the crust is marginal, since most 
of the pressure support is provided by the degenerate elec- 
trons. Still, there are some interesting nuclear-structure ef- 
fects that must be taken into consideration. For example, lack 
of experimental information on the masses of the neutron-rich 


N = 82 isotones generates some uncertainty in the composi- 
tion of the outer crust. Given our limited knowledge of the 
density dependence of the symmetry energy, this affects the 
nuclear composition in the bottom layers of the outer crust 
which, in turn, affects the electron fraction and consequently 
the pressure support. However, the impact of an uncertain nu- 
clear composition on the EOS is very small. Uncertainties in 
the inner crust are more severe because the ground state af- 
ter neutron drip is complex and because no guidance can be 
obtained from experiment. However, whereas the existence 
of these exotic phases influence the transport properties in the 
stellar crust, their impact on the EOS is believed to be rather 
small [54]. Henceforth, we focus exclusively on the study of 
the EOS above the crust-core transition density where the sys- 
tem is uniform and nuclear interactions play a predominant 
role. 


C. Theoretical Constraints 


The saturation of symmetric nuclear matter, namely, the ex- 
istence of a nearly uniform density in the nuclear interior of 
medium to heavy nuclei is one of the hallmarks of the nuclear 
dynamics. By the same token, finite nuclei can only provide 
a very narrow window into the equation of state. Thus, one is 
forced to resort to theoretical approaches to compute the EOS 
below saturation density. Using a Quantum Monte Carlo ap- 
proach, the equation of state of dilute neutron matter has been 
computed with a simplified interaction that matches smoothly 
to the known analytic results at very low densities and pro- 
vides important information at higher densities [167]. The 
choice of interaction has since been refined by adopting a chi- 
ral interaction that is motivated by the underlying symmetries 
of QCD. One of the great virtues of chiral effective field the- 
ory (EFT) is that it provides a systematic and improvable ex- 
pansion in terms of a suitable small parameter [166] [168] [169]; 
for a recent review on the application of chiral EFT to the 
nucleon-nucleon interaction see Ref. [170]. During the last 
decade, enormous progress has been made in our understand- 
ing of the EOS of pure neutron matter by systematically im- 
proving the chiral expansion [171/177]. Although the chiral 
expansion is known to break down above saturation density, 
high-order chiral EFT calculations of pure neutron matter pro- 
vide vital constraints on the EOS at and below saturation den- 


sity [178]. 


D. Laboratory Constraints 


Laboratory experiments on finite nuclei play an essential 
role in constraining the EOS of neutron rich matter in the 
vicinity of saturation density. Although no single technique 
can be used to determine the EOS over the enormous range of 
densities spanned in a neutron star, several methods provide 
insights over a common range of densities, which can then be 
used to test for consistency and for extending the EOS to the 
next density rung. In particular, laboratory experiments pro- 
vide valuable information on the density dependence of the 


symmetry energy. Recall that the symmetry energy is related 
to the energy cost required for converting symmetric nuclear 
matter into pure neutron matter. In particular, the slope of 
the symmetry energy at saturation density is closely related to 
the corresponding pressure of pure neutron matter—a quantity 
that has already been predicted using chiral EFT (171/177). 


The slope of the symmetry energy, a quantity denoted as L, 
may also be extracted from laboratory experiments. Although 
not directly an observable, L is strongly correlated to the thick- 
ness of the neutron skin of heavy nuclei [179/181]. Because 
of the repulsive Coulomb interaction between protons, heavy 
nuclei tend to be neutron rich. For example, the nucleus of 
208Pb contains 126 neutrons but only 82 protons. Hence, the 
root-mean-square radius of the neutron distribution tends to 
be larger than the radius of the corresponding proton distribu- 
tion. The difference between the neutron and proton radii is 
defined as the neutron skin. 


But what controls the thickness of the neutron skin? [182]. 
Although the quantitative value of the neutron skin requires 
detailed calculations, a qualitative explanation of its size and 
its strong connection to the slope of the symmetry energy are 
largely model independent. Surface tension favors the for- 
mation of a spherical liquid drop containing all 208 nucleons. 
However, the symmetry energy increases monotonically in the 
density region of relevance to atomic nuclei. So the symme- 
try energy favors moving some of the excess neutrons to the 
surface where the symmetry energy is lower than in the inte- 
rior. This, however, increases the surface tension. Ultimately 
then, the thickness of the neutron skin emerges from a compe- 
tition between the surface tension and the difference between 
the value of the symmetry energy in the interior relative to 
its value at the surface. Such a difference is encoded in the 
slope of the symmetry energy L. If this difference is large, 
namely, there is a significant energy gain by moving some 
excess neutrons to the surface, then the neutron skin will be 
thick. Conversely, if the energy cost is low, then the excess 
neutrons will remain in the core, leading to a thin neutron skin. 
This suggests a powerful correlation: the larger the value of 
L the thicker the neutron skin [183]. Such a robust correla- 
tion has been verified using a large set of models that have 
been calibrated to reproduce well-known properties of finite 


nuclei [181]. 


Motivated by the possible impact of a terrestrial experi- 
ment in constraining the EOS, a proposal was submitted to the 
Thomas Jefferson National Accelerator Facility more than 20 
years ago. The proposed goal of the experiment was to mea- 
sure the neutron skin thickness of 7°8Pb using parity violating 
electron scattering; denoted by PVES in Fig] By measuring 
the neutron skin thickness of 7°°Pb, one would be able to pro- 
vide valuable information on L. Unlike experiments involv- 
ing strongly-interacting probes, such as pions and protons, 
PVES is clean because it relies exclusively on electroweak 
probes. Indeed, the parity-violating asymmetry emerges from 
the interference of Feynman diagrams involving only the ex- 
change of photons and Z° bosons. Moreover, unlike the pho- 
ton, which is largely insensitive to the neutron distribution, the 
neutral weak vector boson Z° couples strongly to the neutrons, 
making it an ideal probe to map the neutron density. 


The first Lead Radius Experiment (“PREX’’) established 
the existence of a neutron-rich skin in ?°8Pb [184] [185], sug- 
gesting that the pressure of pure neutron matter is “stiff”, 
namely, that the pressure increases rapidly with increasing 
density in the vicinity of saturation density. Unfortunately, the 
impact of the experiment was hindered by unanticipated large 
error bars. A second effort aimed at improving the precision 
of the original experiment was successfully completed and re- 
ported a neutron skin thickness that confirmed—with much 
better precision—that the neutron skin of 2°8Pb is thick 
and, in turn, that the symmetry energy is stiff [126]. Relying 
on the strong correlation between the neutron skin thickness 
of 7°8Pb and the slope of the symmetry energy [181], one was 
able to infer a value for L that systematically overestimates ex- 
isting limits based on both theoretical approaches and exper- 
imental measurements involving strongly-interacting probes. 
Hence, PREX provides the strongest evidence to date in fa- 
vor of a stiff EOS around saturation density. Moreover, be- 
cause neutron-star radii are mostly sensitive to the pressure 
around twice saturation density [186], PREX constrains both, 
the thickness of the neutron skin and, although to a lesser 
extent, the radius of a neutron star [187]. This is a remark- 
able connection involving objects that differ in size by more 
than 18 orders of magnitude [188]. Other laboratory experi- 
ments, such as those that excite the isovector giant dipole res- 
onance or heavy-ion collisions (IVGDR and HIC in Fig/4) are 
also sensitive to the density dependence of the symmetry en- 
ergy. A discussion of the advantages and disadvantages of 
these and many other experimental techniques may be found 


in Refs. [1894191]. 


E. Observational Constraints I 


The first rung in the ladder relies on chiral EFT to determine 
the EOS at low densities. Piggybacking on this rung, labora- 
tory experiments extend our knowledge of the EOS perhaps as 
far as twice saturation density. In turn, measurements of stel- 
lar radii by the Neutron Star Interior Composition Explorer 
(NICER) inform the EOS at even higher densities. NICER, 
launched in June 2017 aboard SpaceX’s Falcon 9 rocket and 
deployed to the International Space Station, is part of NASA’s 
first program dedicated to the study of the exotic structure and 
composition of neutron stars. NICER relies on the powerful 
technique of Pulse Profile Modeling to monitor electromag- 
netic emission from the hot spots located on the surface of the 
neutron star [46][192]. Magnetic fields in pulsars are so strong 
and complex that charged particles that are ripped away from 
the star often crash back into the stellar surface creating hot 
spots, regions within the star that glow brighter than the rest 
of the star. As the neutron star spins, the hot spots come in and 
out of view, producing periodic variations in the brightness, or 
a distinct pulse profile, that are recorded by NICER’s power- 
ful instruments. Because of gravitational light-bending—one 
of the most dramatic consequence of the general theory of 
relativity—hot spots emissions are detected by NICER, even 
when the emission is coming from the “back” of the star. The 
more compact the star, the more dramatic the effect. In this 


manner, NICER provides critical information on the EOS by 
constraining the stellar compactness: 
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where v, M, R, and R, are the escape velocity, mass, ra- 
dius, and Schwarzschild radius of the star, respectively. The 
Schwarzschild radius represents the radius at which the neu- 
tron star would become a black hole, or equivalently, when the 
escape velocity from the star would equal the speed of light. 
For example, a one-solar mass object like our Sun would be- 
come a black hole when its radius reaches about 3 km. 

Before the deployment of NICER, no single neutron star 
had both their mass and radius simultaneously determined. 
This changed in 2019 with the determination of the mass and 
radius of the millisecond pulsar PSR J0030+0451. Two dif- 
ferent groups, one based in Europe and the other one 
in the United States [50], have provided independent—and 
fully consistent—estimates. Both groups found that a neu- 
tron star with a mass of about 1.4 solar masses has a radius 
of approximately 13 km [49] [50]. The second target of the 
NICER mission was the millisecond pulsar PSR JO740+6620 
that was selected because, unlike PSR J0030+0451, its mass 
was previously determined using pulsar timing 33]; see 
Section[IIIG| For this massive pulsar with a mass of more than 
twice the mass of the Sun, the extracted neutron star radius 
was about 12.4km [51]|52]. In tandem, these two pioneering 
observations suggest that the equation of state in the relevant 
density domain is relatively stiff. Moreover, these two mea- 
surements appear to validate an earlier conjecture that sug- 
gests that neutron stars have approximately the same radius 
over a wide mass range [41]. 


FE. Gravitational-Wave Constraints 


The historical detection of gravitational waves emitted from 
the binary neutron star merger GW 170817 has opened a brand 
new window into the Universe [1]. Gravitational waves offer 
a treasure trove of new information that complements elec- 
tromagnetic observations of myriad of astrophysical phenom- 
ena. In particular, GW170817 started to provide answers to 
two fundamental questions in nuclear science: (1) What are 
the new states of matter at exceedingly high density and tem- 
perature? and (ii) how were the elements from iron to uranium 
made? 

Just a few hours after the gravitational-wave detection, 
ground- and spaced-based telescopes identified the associated 
electromagnetic transient, or “kilonova’, that is believed to 
have been powered by the radioactive decay of the heavy el- 
ements synthesized in the rapid neutron-capture process 
' 

Concerning the equation of state, GW170817 provides a 
few critical observables, such as the chirp mass, the mass 
ratio, and the tidal deformability. The tidal deformability 
describes the tendency of a neutron star to develop a mass 
quadrupole in response to the tidal field generated by its com- 
panion star [193] {194]. In the linear regime of small tidal dis- 


turbances, the constant of proportionality connecting the tidal 
forces to the mass quadrupole is the dimensionless tidal de- 
formability A defined as 
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where ky is the second Love number 195) [196]. Although k» 
is sensitive to the underlying EOS, most of the sensitivity is 
contained in the compactness parameter [3] [197/203]. Indeed, 
for a given mass, the tidal deformability scales approximately 
as the fifth power of the stellar radius, a quantity that has been 
notoriously difficult to constrain before NICER and LIGO- 
Virgo [48]. In particular, an estimate by the LIGO- 
Virgo collaboration of the tidal deformability of a 1.4Mgun 
neutron star yields the relatively small value of Ay,4 < 580, 
suggesting that neutron stars are compact objects that are dif- 
ficult to tidally deform [5]. Such a small value for the tidal 
deformability translates into a stellar radius of approximately 
R,4=11.4km [204]. Although consistent within error bars, 
such a low central value for R,.4 differs significantly from the 
R,.4=13km radius extracted from the first NICER observa- 
tion. To our knowledge, GW170817 provides one of the very 
few indications that the EOS is soft. 

Besided the tidal deformability, other observables of rele- 
vance to the EOS are encoded in the gravitational waveform. 
To leading order in the post-Newtonian expansion—an ex- 
pansion in terms of the ratio of the orbital velocity to the 
speed of light—the gravitational-wave signal is determined 
by the “chirp” mass, which involves a linear combination 
of the individual stellar masses [205]. To determine the in- 
dividual masses one must invoke higher orders in the post- 
Newtonian expansion [206]. Although the individual masses 
of GW170817 could not be determined very precisely, the ex- 
tracted upper limit on the most massive of the two stars 
did not challenge the present maximum-mass limit of about 
2.1M.. This situation could have changed dramatically with 
the detection of gravitational waves from the coalescence of 
a binary system GW190814 with the most extreme mass ra- 
tio ever observed: a 23 solar mass black hole and a 2.6 so- 
lar mass compact object [207]. The LIGO-Virgo collaboration 
suggested that GW 1908 14 is unlikely to have originated from 
a neutron star-black hole coalescence, yet left open the pos- 
sibility that improved knowledge of the equation of state or 
future observations could alter this assessment. 

Within a theoretical framework uniquely suited to inves- 
tigate both the properties of finite nuclei and neutron stars, 
it was concluded that the low deformability demanded by 
GW170817 combined with heavy-ion data, make it highly 
unlikely that neutron stars can have masses as large as 
2.6Mz [208]. The absence of very massive neutron stars is 
also consistent with the analysis by Margalit and Metzger who 
argue against their formation based on the lack of evidence of 
a large amount of rotational energy in the ejecta during the 
spin-down phase of GW170817 [209]. Interestingly, the sug- 
gested upper limit of Mmax <2.17Mo is in full agree- 
ment with the recent observation by Cromartie and collabo- 
rators of the most massive neutron star to date 33]. Yet 
the mystery remains unsolved, as the existence of black holes 
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with a mass as small as 2.6M. is also problematic. It has 
even been suggested that GW 1908 14 may have been a second- 
generation merger from a hierarchical triple system: a binary 
neutron star coalescence leading to the 2.6 Ms compact ob- 
ject that ultimately merges with the 23 Mj black hole {210}. 
Finally, the detection of massive neutrons stars from the bi- 
nary coalescence complements pulsar timing techniques that 
provide the most precise information to date on massive neu- 
tron stars and forms the last observational rung in the EOS 
ladder. 


G. Observational Constraints II 


According to Newton’s law of universal gravitation, all that 
can be determined from the motion of two orbiting objects is 
their combined mass. This fact is encapsulated in Kepler’s 
third law of planetary motion: 
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(5) 
where P is the orbital period, a is the length of the semi- 
major axis of the ellipse, and M; and Mp are the individual 
masses of the two orbiting bodies. To determine the individ- 
ual masses one must go beyond Newtonian mechanics and 
invoke general relativity. The most successful observational 
technique used to determine the individual masses is pulsar 
timing. The success of pulsar timing in extracting orbital pa- 
rameters stems from the remarkable precision of pulsars in 
keeping time. These celestial timekeepers often rival the pre- 
cision of the most sophisticated atomic clocks. Millisecond 
pulsars form a distinct group of recycled neutron stars that 
rotate extremely fast. It is believed that the rapid rotation in- 
volves a spin-up during mass accretion from a companion star. 
Further, this class of pulsars generate weak magnetic fields of 
“only” about 100 million Gauss, leading to very small energy 
losses due to magnetic dipole radiation. The fast rotation and 
the weak magnetic fields make pulsars extremely stable clocks 
with minuscule spin-down rates of about 1 second every tril- 
lion years! Long-time observations of millisecond pulsars that 
account for each rotation provide a precise determination of 
the orbital parameters. 

Shapiro delay [211], regarded as the fourth test of general 
relativity, is a powerful technique that has been used in com- 
bination with pulsar timing to infer the mass of the compan- 
ion star, thereby breaking the mass degeneracy encountered 
in Newtonian gravity. The main concept behind the Shapiro 
delay is that even massless particles, such as photons, bend 
as they pass near massive objects because of the curvature of 
space-time. This causes a time delay in the arrival of the elec- 
tromagnetic radiation emitted by the neutron star as its passes 
near the companion star on its way to the observer. In ac- 
counting for every rotation of the neutron star over long peri- 
ods of time, pulsar timing together with Shapiro delay pro- 
vide a highly precise value for the mass of the companion 
star. Once the mass degeneracy has been broken, Kepler’s 
law may then be used to extract the mass of the neutron star. 
This powerful technique was used back in 2010 to determine 


both the mass of the pulsar J1614-2230 and of its compan- 
ion white-dwarf star [30]. In particular, the mass of the neu- 
tron star was found to be close to two solar masses, providing 
the first serious challenge to models that predict soft equa- 
tions of state due to the emergence of exotic states of mat- 
ter. Recently, Shapiro delay was used to measure the mass of 
the millisecond pulsar JO740+6620, that with a mass of about 
M=2.1Mo is the most massive—well measured—neutron 
start to date [33]; although see [34]. Note that the mil- 
lisecond pulsar JO740+6620 is one of the two neutron stars 
that was targeted by the NICER mission and for which a ra- 
dius measurement now exists [51] [52]. 


H. Conclusions 


So what can we conclude from the wealth of information 
that we have gathered during the last five years? Figure[5] 
encapsulates some of the information obtained from labo- 
ratory experiments (PREX-2), electromagnetic observations 
(NICER), and gravitational-wave detections (LIGO-Virgo). 
Predictions are displayed from a set of models that repro- 
duce properties of finite nuclei and generate an equation of 
state stiff enough to support neutron stars with masses of 
at least 2M.. At the lo level, PREX-2 disfavors models 
that predict a neutron skin thickness in *°8Pb of less than 
Rgkin =0.21 fm [125], suggesting that the EOS in the vicinity 
of saturation density is stiff [126]. In turn, the two magenta 
horizontal lines indicate the radius of aM =1.4M. neutron 
star as reported by the two independent analyses of the NICER 
data [49] 50]. Although differing by 18 orders of magnitude in 
size and probing different regions of the EOS, both PREX-2 
and NICER favor a relatively stiff equation of state. However, 
some tension emerges as we incorporate the tidal deformabil- 
ity of aM=1.4M, neutron star extracted from GW170817. 
Combining two analyses [[I] [5], the recommended upper limit 
for the tidal deformability was reported to be 580 at the 90% 
confidence level. As it stands, all models displayed in the fig- 
ure will be inconsistent with such an upper limit. Such a soft- 
ening of the EOS at intermediate densities—bracketed by a 
stiff EOS at both low densities and high densities as suggested 
by PREX-2 and massive neutron stars, respectively, may be 
indicative of a phase transition in the stellar interior [3]. While 
enormously exciting, one must wait for the next generation 
of experimental and observational facilities to validate such a 
scenario. 

Finally, Fig{6] displays the holy grail of neutron-star struc- 
ture, namely, the mass-radius (MR) relation informed by the 
many pioneering developments in the field. Mass constraints 
obtained from electromagnetic observations of the three most 
massive neutron stars to date are shown with the horizontal 
bars [30133]. In turn, the two horizontal arrows in the figure 
indicate constraints on stellar radii obtained exclusively from 
GW 170817, suggesting that the equation of state in the rele- 
vant density region is soft [2)[3]. Such a scenario excludes sev- 
eral models that predict larger stellar radii. By measuring rel- 
atively large radii for PSR JO030 and PSR JO740, NICER has 
relaxed some of the tight constraint imposed by GW170817. 
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Figure 5. Experimental and Observational constraints on theoretical 
models of the equation of state, with its predictions denoted by the 
blue dots. Results are shown for the neutron skin thickness of 208 pp, 
together with the tidal deformability and radius of aM=1.4Mz neu- 
tron star. The combined PREX-II result together with NICER con- 
straints on the stellar radius is depicted by the small (blue) window 
of models allowed [126]. 


NICER results are denoted in the figure by the horizontal error 
bars located at about M=1.4M. and M=2M.Q. The excluded 
causality region was adopted from Ref.[186]. Relative to the 
status of the field five years ago, the progress in tightening the 
MR relation—and consequently the EOS—is simply remark- 
able. 

In summary, due a confluence of pioneering developments 
in a variety of interconnected disciplines, we have entered 
the golden era of neutron stars [213]. Remarkable ad- 
vances are anticipated across all these disciplines during the 
next few years. Among these, the operation of brand new ter- 
restrial facilities that will measure the properties of atomic nu- 
clei at the limits of existence, new telescopes operating across 
the entire electromagnetic spectrum that will determine with 
unprecedented precision the mass-radius relation, and future 
gravitational-wave observatories that will make countless de- 
tections of gravitational waves from compact binary coales- 
cence. Such remarkable developments, combined with ad- 
vances in our theoretical understanding of dense, neutron- 
rich matter, will ultimately lead to a detailed understanding 
of the structure, dynamics, and composition of neutron stars, 
arguably, among the most fascinating objects in the universe. 
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